Manufactured solutions for Elasticity

In [1]:

    

from __future__ import division
from sympy import *
x, y, z, t = symbols('x y z t')
f, g, h = symbols('f g h', cls=Function)
init_printing()


    

In [2]:

    

L = symbols('L')
a1, a2, a3, b1, b2, b3, c1, c2, c3 = symbols('a1 a2 a3 b1 b2 b3 c1 c2 c3')
u0, ux, uy, uz, v0, vx, vy, vz, w0, wx, wy, wz = symbols('u0 u_x u_y u_z\
                                                         v0 v_x v_y v_z\
                                                         w0 w_x w_y w_z')
lamda, mu, rho = symbols('lamda mu rho')


    

In [3]:

    

u = u0 + ux*sin(a1*pi*x/L) + uy*sin(a2*pi*y/L) + uz*sin(a3*pi*z/L)
v = v0 + vx*sin(b1*pi*x/L) + vy*sin(b2*pi*y/L) + vz*sin(b3*pi*z/L)
w = w0 + wx*sin(c1*pi*x/L) + wy*sin(c2*pi*y/L) + wz*sin(c3*pi*z/L)


    

In [4]:

    

def laplacian(U, X=[x, y, z]):
    lap_U = Matrix([sum([diff(U[k], X[j], 2) for j in range(3)]) for k in range(3)])
    return lap_U  


def div(U, X=[x, y, z]):
    return sum([diff(U[k], X[k]) for k in range(3)])


def grad(f, X=[x, y, z]):
    return Matrix([diff(f, X[k]) for k in range(3)])


def rot(U, X=[x, y, z]):
    return Matrix([diff(U[2], X[1]) - diff(U[1], X[2]),
                   diff(U[0], X[2]) - diff(U[2], X[0]),
                   diff(U[1], X[0]) - diff(U[0], X[1])])

def navier(U, X=[x,y,z], lamda=lamda, mu=mu):
    return mu*laplacian(U, X) + (lamda + mu)*grad(div(U, X), X)


def dt(U, t=t):
    return Matrix([diff(U[k], t) for k in range(3)])


    

In [5]:

    

U = [u,v,w]
F = rho*dt(U) - navier(U)

simplify(F)


    

    
Out[5]:





$$\left[\begin{matrix}\frac{\pi^{2}}{L^{2}} \left(a_{1}^{2} u_{x} \left(\lambda + \mu\right) \sin{\left (\frac{\pi a_{1}}{L} x \right )} + \mu \left(a_{1}^{2} u_{x} \sin{\left (\frac{\pi a_{1}}{L} x \right )} + a_{2}^{2} u_{y} \sin{\left (\frac{\pi a_{2}}{L} y \right )} + a_{3}^{2} u_{z} \sin{\left (\frac{\pi a_{3}}{L} z \right )}\right)\right)\\\frac{\pi^{2}}{L^{2}} \left(b_{2}^{2} v_{y} \left(\lambda + \mu\right) \sin{\left (\frac{\pi b_{2}}{L} y \right )} + \mu \left(b_{1}^{2} v_{x} \sin{\left (\frac{\pi b_{1}}{L} x \right )} + b_{2}^{2} v_{y} \sin{\left (\frac{\pi b_{2}}{L} y \right )} + b_{3}^{2} v_{z} \sin{\left (\frac{\pi b_{3}}{L} z \right )}\right)\right)\\\frac{\pi^{2}}{L^{2}} \left(c_{3}^{2} w_{z} \left(\lambda + \mu\right) \sin{\left (\frac{\pi c_{3}}{L} z \right )} + \mu \left(c_{1}^{2} w_{x} \sin{\left (\frac{\pi c_{1}}{L} x \right )} + c_{2}^{2} w_{y} \sin{\left (\frac{\pi c_{2}}{L} y \right )} + c_{3}^{2} w_{z} \sin{\left (\frac{\pi c_{3}}{L} z \right )}\right)\right)\end{matrix}\right]$$

In [6]:

    

L = symbols('L')
h0, hx, hy = symbols('h0 hx hy')
lamda, mu, rho = symbols('lamda mu rho')
alpha, mu_c, xi, gamma, J = symbols('alpha mu_c xi gamma J')


    

In [7]:

    

u = u0 + ux*sin(a1*pi*x/L) + uy*sin(a2*pi*y/L)
v = v0 + vx*sin(b1*pi*x/L) + vy*sin(b2*pi*y/L)
w = 0
h1 = 0
h2 = 0
h3 = h0 + hx*sin(c1*pi*x/L) + hy*sin(c2*pi*y/L)


    

In [8]:

    

def coss(U, H, X=[x,y,z], lamda=lamda, mu=mu, mu_c=mu_c, xi=xi, gamma=gamma):
    f_u = (lamda + 2*mu)*grad(div(U)) - (mu + mu_c)*rot(rot(U)) + 2*mu_c*rot(H)
    f_h = (alpha + 2*gamma + 2*xi)*grad(div(H)) - 2*xi*rot(rot(H)) + 2*mu_c*rot(U) - 4*mu_c*Matrix(H)
    return f_u, f_h


    

In [9]:

    

U = [u, v, w]
H = [h1, h2, h3]
f_u, f_h = coss(U, H)


    

In [12]:

    

simplify(f_u)


    

    
Out[12]:





$$\left[\begin{matrix}\frac{\pi}{L^{2}} \left(2 L c_{2} hy \mu_{c} \cos{\left (\frac{\pi c_{2}}{L} y \right )} - \pi \left(a_{1}^{2} u_{x} \left(\lambda + 2 \mu\right) \sin{\left (\frac{\pi a_{1}}{L} x \right )} + a_{2}^{2} u_{y} \left(\mu + \mu_{c}\right) \sin{\left (\frac{\pi a_{2}}{L} y \right )}\right)\right)\\- \frac{\pi}{L^{2}} \left(2 L c_{1} hx \mu_{c} \cos{\left (\frac{\pi c_{1}}{L} x \right )} + \pi \left(b_{1}^{2} v_{x} \left(\mu + \mu_{c}\right) \sin{\left (\frac{\pi b_{1}}{L} x \right )} + b_{2}^{2} v_{y} \left(\lambda + 2 \mu\right) \sin{\left (\frac{\pi b_{2}}{L} y \right )}\right)\right)\\0\end{matrix}\right]$$

In [13]:

    

simplify(f_h)


    

    
Out[13]:





$$\left[\begin{matrix}0\\0\\- \frac{1}{L^{2}} \left(4 L^{2} \mu_{c} \left(h_{0} + hx \sin{\left (\frac{\pi c_{1}}{L} x \right )} + hy \sin{\left (\frac{\pi c_{2}}{L} y \right )}\right) + 2 \pi L \mu_{c} \left(a_{2} u_{y} \cos{\left (\frac{\pi a_{2}}{L} y \right )} - b_{1} v_{x} \cos{\left (\frac{\pi b_{1}}{L} x \right )}\right) + 2 \pi^{2} \xi \left(c_{1}^{2} hx \sin{\left (\frac{\pi c_{1}}{L} x \right )} + c_{2}^{2} hy \sin{\left (\frac{\pi c_{2}}{L} y \right )}\right)\right)\end{matrix}\right]$$

In [ ]:

    




    

In [ ]:

    




    

In [12]:

    

from IPython.core.display import HTML
def css_styling():
    styles = open('./styles/custom_barba.css', 'r').read()
    return HTML(styles)
css_styling()


    

    
Out[12]:

In [ ]: